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Abstract 

We study the next-to-next-to-leading order (NNLO) evolution of flavour singlet parton 
densities and structure functions in massless perturbative QCD. Present information on 
the corresponding three-loop splitting functions is used to derive parametrizations of 
these quantities, including Bjorken-x dependent estimates of their residual uncertainties. 
Compact expressions are also provided for the exactly known, but in part rather lengthy 
two-loop singlet coefficient functions. The size of the NNLO corrections and their effect 
on the stability under variations of the renormalization and mass-factorizations scales are 
investigated. Except for rather low values of the scales, the residual uncertainty of the 
three-loop splitting functions does not lead to relevant effects for x > 10~ 3 . Inclusion of the 
NNLO contributions considerably reduces the theoretical uncertainty of determinations 
of the quark and gluon densities from deep-inelastic structure functions. 
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1 Introduction 



Structure functions in deep-inelastic scattering (DIS) form the backbone of our knowledge 
of the proton's parton densities — which are indispensable for analyses of hard scattering 
processes at proton-(anti-)proton colliders like Tevatron and the future LHC — and are 
among the quantities best suited for measuring the strong coupling constant a s . The full 
realization of this potential, however, requires transcending the standard next-to-leading 
order (NLO) approximation of perturbative QCD summarized in ref. In fact, the 
next-to-next-to-leading order (NNLO) subprocess cross-sections ('coefficient functions') 
for DIS have been calculated some time ago [0, [| (see also refs. [Q for the related Drell- 
Yan process), but the corresponding three- loop splitting functions governing the NNLO 
scale dependence ('evolution') of the parton densities have not been completed so far (see 
ref. || for an up-to-date progress report). 

In a previous paper || we have studied the flavour non-singlet sector of the DIS 
structure functions. It turned out that the available incomplete information on the non- 
singlet splitting functions ]7|, [|, || [HJ facilitates the derivation of approximate expressions 
which are sufficiently accurate over a rather wide region of parton momenta and the 
Bjorken variable x. In the present paper, we complement those results by corresponding 
effective parametrizations for the flavour-singlet sector based on the partial information 

-3 



of refs. P, FT, 12, 13], thus paving the way for promoting, even though only at a;>10 



global analyses |L4|, |T^, [16| of DIS and related processes [|J to NNLO accuracy. 



This paper is built up as follows: In Sect. 2 we set up our notations by recalling the 
general framework for the evolution of singlet parton densities and structure functions in 
massless perturbative QCD. The expansions of the corresponding splitting functions and 
coefficient functions are written down up to order for arbitrary values of the renormali- 
zation and mass-factorization scales. In Sect. 3 we present compact, but very accurate 
parametrizations of the exactly known fl7| , [| |3|, but partly rather lengthy expressions for 



the two-loop singlet coefficient functions. In Sect. 4 the available constraints |8|, PI |13 
on the three-loop singlet splitting functions are utilized for constructing approximate 
expressions for the x-dependence of these functions, including quantitative estimates of the 
remaining uncertainties. In Sect. 5 we assemble all these results and quantify the impact of 
the NNLO contributions on the evolution of the singlet quark and gluon densities and on 
the most important singlet structure function, F 2 . We address the range of applicability 
of the present approximate results and the improvement of the theoretical accuracy of 
determinations of the parton densities at NNLO. Finally our results are summarized in 
Sect. 6. Mellin-iV space expressions for our parametrizations of the two- loop coefficient 
functions of Sect. 3 can be found in the appendix. 
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2 General Framework 



We start by outlining the general formalism for the NNLO evolution of flavour- singlet 
parton densities and structure functions. The singlet quark density of a hadron is given by 



»=i 



Here qi(x, fij, fj%) and q~i(x, fj/j, /jl 2 .) represent the number distributions of quarks and anti- 
quarks, respectively, in the fractional hadron momentum x. The corresponding gluon 
distribution is denoted by g(x, /i 2 ). The subscript i indicates the flavour of the (anti-) 
quarks, and Nf stands for the number of effectively massless flavours. Finally fi r and 
fif represent the renormalization and mass-factorization scales, respectively. The singlet 
quark density (2.1) and the gluon density are constrained by the energy- momentum sum 
rule 

/•i 

E(x,fj? f ,t£)+g(x,ti},t£))\ = 1 . (2.2) 



i 

dx x 



o 



The scale dependence of the singlet parton densities is given by the evolution equations 

d\nn) ciln/4 \g ) \ V gq V gg J \g ) 
where ® stands for the Mellin convolution in the momentum variable, 

[a®b\{x] = [ ld ^-a{y)b(-) . (2.4) 
Jx y \yj 

As in some other equations below, the dependence on x, [if and \i r has been suppressed 
in Eq. (2.3). The splitting function V qq can be expressed as 

V gq = V+ s + V PS (2.5) 

with 

^ps = N f {V s qq + V s qq ) . (2.6) 

Here V^ s is the non-singlet splitting function discussed up to NNLO in ref. 0, and the 
0{a 2 s ) quantities V qq and V^ q are the flavour independent ('sea') contributions to the 
quark-quark and quark- ant iquark splitting functions V qiqk and V qiqk , respectively. The 
non-singlet contribution dominates Eq. (2.5) at large x, where the pure singlet term P PS 
is very small. At small x, on the other hand, the latter contribution takes over as xV PS , 
unlike xV^ s , does not vanish for x— >0. The gluon-quark and quark-gluon entries in (2.3) 
are given by 

P» = NfV qi9 , V gq = V m (2.7) 
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in terms of the flavour- independent splitting functions V qi9 = V qi9 and V gqi = V gqi . With 
the exception of the lowest-order approximation to V qg , neither of the quantities xV qg , 
xV gq and xV gg vanishes for x— >0. 

The NNLO expansion of the splitting-function matrix P in Eq. (2.3) readsQ 

P(x,a.{t4),L R ) = a s P (0 \x) 

+ a 2 s (pM(x) - f3 L R P {0 \x)) (2.8) 
+ a 3 s (P®(x) - 2p L R P^\x) - {p x L R - (3 2 L 2 R }P^°\x)) +... . 

In Eq. (2.8) and in what follows we use the abbreviations 



An 

for the running coupling, and 



(2.9) 



L M = \n% , L R = ln^ (2.10) 

for the scale logarithms. The one- and two-loop matrices P^°\x) and P^\x) in Eq. (2.8) 
are known for a long time 0, see also ref. [18| for the solution of an earlier problem in 



the covariant-gauge calculation of PjQ ■ The three- loop quantities P^\x) are the subject 
of Sect. 4. The consistency of the evolution equations (2.3) with the momentum sum rule 
(2.2) imposes the following constraints on the second moments of P^\x) : 

Pi l q ] (N) + P^(N) = 

P$(N) + P$(N) = for iV = 2, (2.11) 

where 

a(N) = [ dxx N ~ 1 a(x) . (2.12) 
Jo 

The constants j3i in Eq. (2.8) represent the perturbative coefficients of the QCD /9-function 

0(a.) = -E«i +2 A • (2-13) 



dln/i^ 

The coefficients /3 , (3\ and @ 2 required for NNLO calculations can be found in refs. M 
and [jHJ, respectively. Finally the L R terms in Eq. (2.8) are obtained from the expression 
for fM r — /j,f by inserting the expansion of a s (A*/) in terms of a s (/x^), 

a s (fi}) = a s (/i r 2 ) - (3 Q L R a 2 s {£) - [f3 x L R - f3 2 L 2 R } a 3 s (tf.) + ... . (2.14) 



1 Notice that the convention adopted for P^> here and in ref. differs from ref. Q by a factor of 2 l+1 
due to the choice of a s = a s /(47r) instead of a s /(27r) as the expansion parameter in Eq. (2.8), and by a 
factor 1/2 from ref. [FJ due to the choice of fj% instead of fj,f for the differentiation in Eq. (2.3). 
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The singlet structure functions F aj s, a = 1, 2, are in Bjorken-x space obtained by the 
convolution (2.4) of the solution of Eq. (2.3) with the corresponding coefficient functions, 



VaF a ,s(x, Q 2 ) = [C atq (a 8 , L M , L R ) ® E(/xJ, fi 2 ) + C a ^ g (a s , L M , L R ) <g> g(jj? f , // 2 )] (x) 

= C a (a 8 ,L M ,L R ) ® q(ii 2 f ,i4) (x) . (2.15) 

Here the electroweak charge factor is included in r] a , e.g., r/i = 2 (e 2 ) -1 and r^ 2 = (x (e 2 )) -1 
for electromagnetic scattering, with an average squared charge (e 2 ) = 5/18 for an even 
iVy. Up to third order in a s the expansion of the coefficient functions takes the form 

C a (x,a s (/4.),L M ,L R ) = C^\x) + a s C^\x,L M ) + 

a 2 C^(x,L M ,L R ) + a 3 s C^(x,L Ml L R ) + ... (2.16) 

' cf(x) + Eal (<£«>(*) + £ c^(x)L^) + . . . , 

i=l V m=l ) 

where c£'(x) = ( c^(x) , c^(x) ) = ( 6(1— x) , ) represents the parton-model result. The 
first-order corrections cft(x) can be found in ref. M; the two-loop coefficient functions 
c^ q (x) and c£> g (x) computed in refs. [!| are briefly addressed in Sect. 3. The terms up to 
order a 2 in Eq. (2.16) contribute to F a in the NNLO approximation. The a 3 terms only 
enter for dF a /d\nQ 2 , hence the scale-independent three-loop quantities c^\x) are not 
needed here. The coefficients c^' m '(x) in the second part of Eq. (2.16) can be conveniently 
written in a recursive manner as 



r (l,l) 










r (2,l) 


= ci°) ® P« + c« 


®(p(°)-A)i) 






r (2,2) 


= IcW) ® (p(°) - 
2 a v 


A)i) 






r (3,l) 


= 4°) ® p( 2 ) + c « 


® (P (1) -fti) 


+ 4 2) < 


g, (p(0) _ 


-(3,2) 
O o 


= ^{ei 1 ' 1) ®(P (1) - 


- ftl) + ci 2 - 1 ) ® 


(p(0) 


-2A)1)} 


r (3,3) 


= ® ( P (o) _ 


2/3 l) • 







(2.17) 

Eqs. (2.17) can be derived by identifying the results of the following two calculations 
of dF a> s/dlnQ 2 at /i 2 = /i 2 = Q 2 : (a) with the scales identified in the beginning (using 
Eqs. (2.3) and (2.13), and (b) with the scales set equal only at the end (after differentiating 
the logarithms in Eq. (2.16) ). Finally the coefficients C^ 2 ' and for /if ^ fi r in 
Eq. (2.16) are obtained from these results by Eq. (2.14): 

C^(x, L M , L R ) = C£\x, L M , 0) - [3 L r c£\x, L M ) 

C®(x,L M ,L R ) = C^(x,L M ,0)-2(3 L r C^>(x,L M ,0) (2.18) 

-{p 1 L R -^L 2 R }C^\x,L M ) . 
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The above-mentioned calculations of fa [0 an d (x) |2], [| have been performed in 
the MS renormalization and factorization schemes. The other factorization scheme widely 
used in NLO analyses of unpolarized deep-inelastic scattering is the so-called DIS scheme. 
In this scheme the expression for F 2j s is reduced, for the choice /J, 2 . — fj, 2 — Q 2 adopted 
for the rest of this section, to its parton-model form 772 -F^s^, Q 2 ) = S(a;, Q 2 ), i.e., 

4%) = c®{x) = for I > 1 . (2.19) 

Here and in the following DIS-scheme quantities are marked by a tilde. The corresponding 
singlet parton densities can be defined by 



Z® q ■ (2.20) 



V 9 J V ~ C2 '" " C2 .f / V 9 J 
The upper row of the transformation matrix Z is fixed by Eq. (2.19). The lower row is 
constrained only by the momentum sum rule (2.2) which fixes its form for the second 
moment, N = 2. The extension of this constraint to all N, 

Z(N)+g(N) = E(N)+g(N) , (2.21) 



in Eq. (2.20) is the obvious generalization of the standard NLO convention of ref. |2(J 
to all orders. The DIS-scheme splitting functions for ft 2 = ft 2 = Q 2 can be expressed in 
terms of their MS counterparts by differentiating of Eq. (2.20) with respect to Q 2 : 

P = (z ®P + P^j® Z" 1 . (2.22) 

Insertion of the expansions (2.8) for Lr = 0, (2.13) and 



Z(x,a s ) = 1 + = 1 + E«i ( <? (0 I (2 ' 23) 



then yields 

P = a s P 



+ a 2 s (P« + [Z« P<°>] - /? .Z (1) ) (2.24) 
+ (p( 2 ) + [Z< 2 >, P (0) ] + [Z«, P (1) ] - [Z«,P (0) ] ® ^ ( 
+ faZ^ g, Z « - 2/3 Z( 2 ) - &Z«) + . . . , 



where 

[Z w , P (m) ] ee Z w <g> P (m) - P (m) ® Z (0 . (2.25) 

The coefficient functions and splitting functions for fi 2 ^ fi 2 ^ Q 2 can by obtained from 
Eqs. (2.19) and (2.24) using Eqs. (2.16)-(2.18) and (2.8), respectively. In the present 
article Eq. (2.24) will be only employed in Sect. 4 for transforming a small- a; constraint 
obtained for P^ 2 ' into the MS scheme adopted for the rest of this paper. 
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3 The 2-loop singlet coefficient functions 



Beyond the first order, 1 = 1, the exact expressions for the scale independent parts c$(x) 
of the coefficient functions (2.16) are rather lengthy and involve higher transcendental 
functions.^ Moreover, the convolutions entering the expressions (2.17) for fi 2 ^ Q 2 (given 
explicitly for I = 2 in ref. JJ) and the factorization scheme transformation (2.24) become 
increasingly cumbersome at higher orders. The latter complications can be circumvented 
by using the moment-space technique [^(], p^j , as the convolution (2.4) reduces to a simple 
product in iV-space. This technique, however, requires the analytic continuation of all 
ingredients to complex values of N in Eq. (2.12), which itself becomes rather involved 
already for the exact expressions for c^ 2 \x) [23]. Hence it is convenient, for both x-space 



and iV-space applications, to dispose of compact parametrizations of these quantities. 

For the non-singlet components c^ NS (x) [|TJ] and c^ s (x) of 

cg(aO = <&(*) + c%(x) , (3.1) 

where Fl = F 2 — 2xFi, such parametrizations have been given in ref. |J in terms of 

L = lux , Li = ln(l - x) (3.2) 

and simpler functions of x. The pure singlet pieces in Eq. (3.1) can be approximated by 

c g) ps (x) = N f {(15.94 -5.212 x) (l-x) 2 L 1 + (0.421 + 1. 520 x)L 2 

+ 28.09 (1 - x)L - (2.370 aT 1 - 19.27)(1 - x) 3 } , (3.3) 

c?-ps( x ) = Nf { — 0.101 (1 — x)L\ — (24.75 — 13.80 x)Li_Lq + 30.23LiL 

+ 4.310 Ll - 2.086 L 2 + 39.78 L + 5.290(2r a - 1)} . (3.4) 

The corresponding gluon coefficient functions can be parametrized as 

cf g (x) = N f {(94.74 -49.20 x)(l - x)L\ + 864.8(1 - x)L x + 1161 xL x L Q 

+ 60.06^ + 39.66 (1 - x)L - 5.333 (aT 1 - 1)} , (3.5) 

d®( x ) = N f {(6.445 + 209.4 (1 - x))L 3 - 24.00 L\ + (1494 x" 1 - 1483)Lx 
+ LiL (-871.8 U - 724.1L ) + 5.319L 3 - 59.48 L 2 - 284.8 L 
+ 11.90X" 1 + 392.4- 0.285(1 - x)} . (3.6) 

Note the small 5(1 — x) term in the parametrization of c^, which is of course absent 
in the exact expression, but useful for obtaining high-accuracy convolutions here. The 
parametrizations (3.3)-(3.6) deviate from the exact results by no more than a few permille. 
2 A recent recalculation of all these functions j2l) has yielded full agreement with refs. |2[ [1. 
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4 The 3-loop singlet splitting functions 



Only partial results have been obtained so far for the O(a^) terms P^\x) of the singlet 
splitting functions (2.8). As the full x-dependence is known just for the CaNj part of Pf® 
| |13| |, the backbone of the available constraints are the lowest four even-integer moments 
(2.12) of these terms, 

P%\N) for N = 2, 4, 6, 8 , (4.1) 
calculated in ref . M . This is one moment less than determined for the non-singlet quantity 



3 (2)+. 
NS 1 



X. 



J/J K, due to the greater technical complexity of the singlet computation. 



The 1/[1 — x]+ soft-gluon contributions to the gluon-gluon splitting functions P$ and 



Pfjg' are related to their quark-quark counterparts by 



x 



x=l 



Ca 
C F 



x)P<g{x) 



x=l 



C, 

c f 



:i-x)p$i( X ) 



x=l 



(4.2) 



(P PS in Eq. (2.5) vanishes at x = 1). This relation also applies to the leading-Afy terms 
|, |T§ of and P^l We will assume that Eq. (4.2) holds generally for I = 2. 

(2) 

The small- a; behaviour of P^j (x) is not constrained by Eq. (4.1). The leading x— >0 
terms oc Mnx have been determined from small-x resummations, however, for Pjjjp , P^ 



and Pgg- Specifically, the gluon-quark and quark-quark entries read |TT 



3 (2) 
qg,x^0 

qq,x^>0 



896 
~27~ 



C\N f 



lnx 



x 



p (2) 

C A q9 ' x ~*° 



O 

0(1 



(4.3) 
(4.4) 



The corresponding gluon-gluon result can be inferred [[24]] from the larger eigenvalue of 
P x 2 ^Lo completed in ref. in a scheme equivalent to the DIS scheme up to 3-loop order: 



P. 



(2) 



gg,x->0 



6320 - 864 C(3) - 1584 C(2) + N f 



1136 



96C(2) 



In x 



x 



0(1 

X 



, (4.5) 



where ((2) = tt 2 /6 and £(3) ~ 1.202057. For this contribution the scheme transformation 
(2.24) schemes simplifies to 



5(2) 
gg,x^0 



3(2) 
gg,x->0 



p(0) 
99 



,(2) 
-2,9 



z->0 



+ - 



(4.6) 



yielding 



P, 



(2) 



gg,x^0 



6320 - 864 C(3) - 1584 C(2) + AT 



/4720 32 



C(2) 



lnx 



.r 



+ O 



• (4-7) 



Eqs. (4.6) and (4.7) are most easily derived via moment-space, see the appendix. 
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In the following the above information will be used for approximate representations of 

4 2) (x) = P®(x) + N f P®(x) + N 2 f P^(x) , (4.8) 

where the iV/-independent terms are absent in P^ and P^> according to Eqs. (2.6) and 
(2.7). For each of the ten remaining contributions P^ m we employ the ansatz 



P§U X ) = E A nfn(x)+fe(x) . (4.9) 



,(2) 

.w - 

n=l 

Except for the 5(1 — x) terms in Pgg] m , the basis functions f n are build up of powers 
of ln(l — x), x, and In a;. Generally we choose the functions f n to include at least one 
function peaking as x—> 1 (except for Ppg, of course), one rather flat function (a small non- 
negative power of x), and one function peaking as x— >0. For each choice of the building 
blocks /„, the coefficients A n are determined from the four linear equations provided by 
the moments (4.1) of ref. || after taking into account the additional information (4.2)- 
(4.4) and (4.7) collected in the function f e in Eq. (4.9). The residual uncertainty of the 
functions Plj' m is estimated by varying the choice of the basis functions f n . Finally two 
approximate expressions are selected for each of the six quantities P g ^ , Pggo and P^\ 
which are representative for the respective present uncertainty band. The Nj pieces P^ 2 
are smaller in absolute size and uncertainty than those contributions, hence for them it 
suffices to select just one central representative. 

Before proceeding to our approximations of the coefficients in Eq. (4.8), it is 

worthwhile to discuss another constraint on the singlet splitting functions which has 
served as an important check for the calculations of the NLO quantities P^\x) in the 



unpolarized [Q] as well as in the polarized case [25|, |2B[ . This is the relation 



pSvr(x) + P£,nn(x) ~ PgM x ) ~ P Sdr(^) = (4.10) 



for the choice 



C A = N c = C F = N s (4.11) 

of the colour factors leading to a jV=l supersymmetric theory. As indicated, Eq. (4.10) 
holds in the dimensional reduction (DR) scheme respecting the supersymmetry, but not in 
the usual MS scheme based on dimensional regularization. The terms breaking the relation 
(4.10) in the MS scheme for I = 1 are, however, much simpler than the functions P^\x) 
themselves. They involve only In 1 x and powers of x, however including the leading small-a; 
term 1/x |[27|| . It should also be noted that Eq. (4.10) does not in all cases introduce a 
relation between different functions P/j . For instance, the leading large-x terms of Pj$ 
and Pq]\ ln 2 (l — x) and ln(l — x), cancel already on the level in the individual splitting 
functions for the choice (4.11), not only after forming the combination (4.10), as can be 



readily read off from the tables in refs. fl28|, E9 . 
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The transformation of the splitting functions from the DRP] to the MS scheme is 
presently unknown for I = 2 (its derivation would require repeating the calculation of 
the finite terms of the two-loop operator matrix elements of ref. |HJ in the DR scheme). 
It is expected, however, that the terms breaking Eq. (4.10) in the MS scheme at I — 2 
are again comparatively simple and do not include terms like ln 3 (l — x), ln 4 (l — x), and 
In 4 x. Under this assumption that relation does provide some constraints on the functions 

(2) 

Pij ■ In view of their limitations discussed above for the NLO case, however, these few 
additional constraints do not justify moving from Eq. (4.8) to a complete colour-factor 
decomposition involving twenty-three instead of ten unknown functions Pij' m . 

We now illustrate the details of our approximation procedure outlined above for the 
case of the gluon-quark splitting functions P qg dominating the small-x evolution of the 
singlet quark density. While the lowest-order result Pj® contains no logarithms, terms 
up to ln 2 (l — x) and In 2 x occur in PW. Hence we expect contributions up to In (1 — x) 
and In 4 x in the three-loop splitting function Pj^ due to the additional loop or emission. 

(2) 

Thus a reasonable choice of the trial functions for Pq 9 \ is given by 

fi(x) = ln 2 x or ln 3 x 

f2{%) = 1 or x 

fz(x) = ln(l— x) or ln 2 (l — x) 

f 4 (x) = ln 3 (l - x) or ln 4 (l - x) 

/.(«) = -«*C. W/ (!=£ + !) , A = 0or4. (4.12) 

As for all other cases in which the - In x terms are known, we do not include the sublead- 
ing contribution 1/x in the functions f n , but vary its coefficient by hand up to a value 
suggested by the expansion of the moment-space expressions P^\n) and P^\n) around 
N = 1 p2| . Four of the thirty-two combinations resulting from Eq. (4.12) (those combin- 
ing fi = In 2 x with f 2 — 1 and A = 4) lead to an almost singular system for the coefficients 
A n which causes overly large oscillations of the corresponding results for pj® x (x). The re- 
maining twenty-eight approximations, of which four are selected for further consideration, 
are displayed in Fig. 1. 

Our corresponding ansatz for the Nj part reads 

fi (x) = In x or In 2 x or In 3 x 
/ 2 (^) = 1 or x 

fz{x) = ln(l — x) or ln 2 (l — x) or ln 3 (l — x) 

U{x) = - and f e (x) = . (4.13) 
x 



3 Actually the construction of a consistent DR scheme respecting the supersymmetry at 3-loop level 
has not been performed so far [M. Consequently Eq. (4.10) has not been proven up to now for l> 1. 
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Figure 1: Approximations of the Nj part Pqg\ of the three-loop splitting function P^(x) 
in Eq. (2.8), as obtained from the four moments (4.1) by means of Eqs. (4.9) and (4.12). 
The full and dashed curves represent the functions selected for further consideration. The 
upper group of curves in the right plot is for A = in Eq.(4.12), the lower group for A = 4. 



Here and for all other Nj terms the \jx contribution is included in the functions f n , as 
no information about its magnitude is available. Fifteen of the eighteen combinations 
resulting from Eq. (4.13) are shown in the left part of Fig. 2. The three combinations 
involving f\ = In 3 x and f2 = x are rejected for the same reason as the four expressions 
for pf£\ mentioned above. 

In the right part of Fig. 2 the approximations selected for and P^ 2 have been 
combined into Eq. (4.8) for Nf = 4 and convoluted with a typical example for the gluon 
density of the proton. Needless to say that it is this combination P^> ®g (divided by g(x) 
for display purposes in Fig. 2) and not P^\x) itself which enters the evolution equations 
(2.3) and thus the structure functions (2.15). At moderately large x, 0.1 <x<0.7, the 
spread of the convolutions in Fig. 2 is small, typically about 5%, much smaller than that 
for the splitting function itself, as the convolution (2.4) tends to wash out the oscillating 
differences of Fig. 1 to a large extent. 

The large-x behaviour of the splitting functions is very important also at small x, as 
can be inferred from a comparison of Figs. 1 and 2: while P^{x) is dominated at x < 10~ 2 
by the \ lnx and 1/x terms, this does by no means hold for [P^ ®g](x). For instance, the 
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X X 

Figure 2: Left: as Fig. 1, but for the Nj contribution Pgf j2 ( x ) using Eqs. (4.9) and (4.13). 
Right: convolution of the approximations A-D of Fig. 1, combined for Nf = 4 with the 
result selected for P^ ' 2 m the left part, with a typical gluon density of the proton. 

leftmost cross-over of the approximations C (A = 4) and D (A = 0) for P^> is at x = 0.03, 
but the convolution for the function D does not exceed that for the approximation C 
above x = 10~ 3 . An even more striking illustration can be obtained by keeping only the 
small-x function f e in Eq. (4.12) for P^i. The convolution then yields (P^ ®g)/g — 4000 
and —6000 for A = and A = 4 at x = 10~ 3 , respectively, instead of the spread between 
about 700 and 1300 found after taking into account the large-x constraints (4.1). 

Except close to the unavoidable cross-over points, the functions A and B are represen- 
tative for the uncertainty bands of Figs. 1 and 2. Hence these approximations are chosen 
as our final results for the gluon-quark splitting function: 

pj;fHx) = -19.5515 L\ + 707.438 L x + 2300.986 + 814.928 hi - — L ° + 4 

3 x 

P^f(x) = 10.8972 L\ — 315.331 L\ + 902.843 x — 1054.09 L\ — — (4.14) 

O Ob 

with L and L\ as defined in Eq. (3.2). The selected approximation for the Nj piece reads 
Pj® 2 (x) = 3.769 L\ — 59.176 x + 0.244 L + 1.079 - . (4.15) 

Ob 

As for all other cases listed below, the average \[A + B] represents our central result. 
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Figure 3: The large-x behaviour of our selected approximations of the three-loop splitting 
functions pj< (x) for Nf = 4. is obtained by adding the pure singlet (PS) term (also 
shown separately) to the non-singlet contribution of ref. || according to Eq. (2.5). 

We now turn to the other splitting functions PpQ, P^ and Pj® . Here we confine 
ourselves for brevity to the final results corresponding to Eqs. (4.14) and (4.15) and 
some brief remarks on the individual cases. The resulting Nf = 4 approximations are 
graphically displayed in Fig. 3 for large x and in Fig. 4 for moderately small x. 

The approximations selected for the pure singlet splitting function are given by 



Pps, 1 ( x ) 
and 



[1 - x) (-93.265 L x + 357.924 x) + 543.482 L 2 + 9.864 L\ 



3584 1 



27 x 



[L + 4(1- x)} 



(4.16) 



;i - x) (37.395 L\ - 210.424 x) - 171.7 L\ - 48.862 Lq 



3584 L 



o 



27 x 



PpS, 2 ( x ) 



1 - x) (-3.999 Li - 12.541 ar) - 8.852 L + 3.445 -(1 - x) 



x 



(4.17) 



12 



The lowest-order pure-singlet contribution Pp$ vanishes like (1 — x) 3 as x— >1. The one 
additional loop or emission in Ppg will introduce logarithms up to ln 2 (l — x), but keep 
Pps( x = 0) = 0. We have implemented this feature into all the basis functions f n in 
Eq. (4.9), e.g., the large-x logarithms are introduced as (1 — x) ln a (l — x). The detailed 

(2) 

x — > 1 behaviour of Pp S is however irrelevant in view of the large-x dominance of the 
non-singlet contribution obvious from Fig. 3. 

The corresponding results for the gluon-gluon splitting function read 

P^q(x) = 2560 . 1 . + 3870.26 5(1 - x) + 1292.56 L x 
(1 — x) , 

- 14903.16 - 3667.22 Ll + 2675.85 -(L + 4) (4.18) 

x 

P^o(x) = 3031 M 1 3 + 5622.22 5(1 - x) - 14514.35 x 2 



(1-x) 



+ 



Pf^(x) = -427.5 - - 570.4 5(1 - x) + 2529.794 x 2 



+ 643.44 + 13565.99 Ll + 2675.85 — 

x 



1 



1605.009 - 784.828 L% + 157.18 -(L + 4) (4.19) 



Pf^(x) = -405.0 . 1 . - 539.16 5(l-a;) - 929.17 L x 

(1 — x) , 



1345.962 x + 21.917 Ll + 157.18 — 

x 



and 



1 G 1 1 

— + 6.575 5(1 - x) - 5.056 x - 9.904 L + 2.969 - 

9(l — a;), x 



The approximate 1/[1 — x] + coefficient in Eqs. (4.18) and (4.19) have been taken over, 
according to Eq. (4.2), from the results (4.9)-(4.11) for the non-singlet splitting functions 
in ref. . The corresponding contribution to Pg g \ is an exact leading- Nf result of ref . Jl^| , 



(2) 

as the second colour factor of P™ 2 n °t determined there does not contain a 1/[1— x] + term. 
The iVj and Nj approximations have been combined in such a manner as to maximize 
the error band where the present uncertainties are largest, i.e., at small x. 

Finally the following expressions are chosen for the quark-gluon splitting function 

P^o(x) = 3.040 Lf + 1157.76 + 2357.73 L 2 + 291. 76 — 

P gq ] ,o(x) = -6.461 Ll - 1789.06 x 2 + 2260.38 - 140.18 , (4.21) 
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Figure 4: As Fig. 3, but for the small- a; behaviour of P^\x). The difference between Pf® 
and Pp<) is negligible for x < 0.1, hence the latter quantity is not separately shown here. 



P^f(x) = 26.2717 L{ + 148.036 L{ - 549.815 x + 89.769 L 3 + 70 



x 



0.1995 L\ - 37.93 L x - 179.644 - 312.616 L\ + 35 



L 



x 



and 



P® 2 (x) = -2.728 L x - 10.217 x- 3.566 L -4.207 



x 



(4.22) 
(4.23) 

Here the coefficient of leading small- a; term Mnx is unknown. The leading x — > con- 
tribution oc 1/x to the LO splitting function Pf® is related to the corresponding term 
of Pf® by a factor Cf/Ca- The same holds, up to two percent, for the iV/-parts of the 
NLO quantities and P$ (but not for the Nj parts). Therefore we have included 



^ In x in the functions /„ in Eq. (4.9) for P gg ) but varied its coefficient by hand for ^ 
between Cf/(2Ca) and Cf/Ca relative to the corresponding term of P 9 ^\ in Eq. (4.19). 
The resulting larger small-x uncertainty of Pjff is clearly visible in Fig. 4. While the 
uncertainty of the other splitting functions relative to the central results not shown in the 
figure does not exceed ±30-40% at x — 10~ 3 , it amounts to about 100% for P@\ 
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5 Numerical results 



We now proceed to illustrate the numerical effects of the NNLO contributions on the 
evolution of the singlet parton densities fj,V) and g(x, (i 2 ) and on the singlet structure 
function F 2) ,s(:r, Q 2 ). Specifically, we will consider the scale derivatives q = dlnq/dlafj, 2 ^, 
q = S, g, at a fixed reference scale [if = /i/,Q , and F 2i s and its Q 2 -derivative at Q 2 = fi 2 Q . 
Except for the last part of our discussions below, these illustrations will be given for the 
fixed, i.e. order-independent, initial distributions 



xJ2(x, /i 2 ^) 



0.6x _a3 (l -x) 3 ' 5 (l + 5x - 8 ) 
1.0x~ a37 (l-x) 5 . 



(5.1) 



These expressions agree well with standard NLO parametrizations ]TJ|, |TB|, [TjJ, at our 
choice fij ~ 30 GeV 2 . Likewise we employ, for the time being, 



0.2 



(5.2) 



irrespective of the order of the expansion, corresponding to a s (M|) ~ 0.116 beyond LO. 
All results below refer to the MS scheme for Nf = 4 massless quark flavours. See ref. 
for a recent discussion of charm mass effects in F%. 
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Figure 5: The perturbative expansion of the scale derivative, S = dlnH/dln fi 2 , of the 
singlet quark density at fj, 2 = /i 2 ~ 30 GeV 2 . The initial conditions are specified in 
Eqs. (5.1) and (5.2). 
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Figure 6: As Fig. 5, but for the gluon density g(x, /i 2 ). The spikes close to x — 0.1 in 
the right parts of both figures are due to zeros of the NLO predictions. 

The LO, NLO and NNLO approximations to the evolution (2.3) of the singlet quark 
and gluon densities are compared in Fig. 5 and Fig. 6, respectively, for the standard choice 
fi r = fif of the renormalization scale. Here and in the following figures the subscripts A 

(2) 

and B refer to the approximate expressions for the functions P^ (x) derived in the previous 
section; the central results i(NNLOA + NNLOs) are not shown separately. For the input 
(5.1) and (5.2), the singlet quark derivative E is dominated by the P qq <S> £ contribution 
at large x, x > 0.3, and by radiation from gluons at small x, x < 0.03 (here the NLO 
corrections are very large, partly due to the absence of 1/x terms in P^ discussed above 
and below Eq. (2.7) ). The gluon derivative g, on the other hand, is mainly driven by the 
P gg <S> g term, but the radiation from quarks is non-negligible over the full x-range. 

In both cases the NNLO corrections are small at large x, amounting to about 2% 
for E and less than 1% for g. The corresponding NLO contributions are 12 — 20% and 

(2) 

2 — 4%, respectively. The present residual uncertainties of P^ are completely immaterial 
in this region of x. The NNLO effects and their uncertainties increase towards very small 
x-values, reaching about (7.5 ± 2)% for £ and (3 ± 2)% for g at x = 10~ 3 . Recall that 
these numbers refer to /ij ~ 30 GeV 2 . At lower scales the small-x shapes of the quark 
and gluon densities are flatter than in Eq. (5.1). Together with a larger a s this leads to 
larger small-x uncertainties. For example, at jjj ~ 3 GeV 2 (corresponding to a s ~ 0.3) 
they reach about ±6% for x = 10~ 3 and fall below ±2% only for x > 4 ® 10~ 3 . 
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The renormalization scale dependence of the NLO and NNLO predictions for the 
derivatives £ and g is shown in Fig. 7 and Fig. 8, respectively, for six representative 
values of x. In both figures /i r is varied over the rather wide interval | fj!f < fif. < 8 fj/j 
corresponding to 0.29 > a s ([i 2 r ) > 0.15 for the input (5.2). The present approximate 
NNLO predictions prove sufficient for a marked improvement on the NLO results, except 
for the bins x = 0.05 (where the derivatives of both the singlet quark and gluon densities 
are small) and, in the gluon 10 3 . As illustrated above the pjp parametrization 

uncertainties of £ and g, which are enhanced at small jj, r due to the larger values of the 
coupling constant, are rather comparable at small x. The present difference between the 
two cases rather stems from the considerably better scale stability of g at NLO. 
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Figure 7: The /^-dependence of the NLO and NNLO predictions for the singlet quark 
derivative S = d In E/d In /x^ at fi 2 f = ~ 30 GeV 2 for six typical values of x. 
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Figure 8: As Fig. 7, but for the gluon derivative g = dlng/dlnfij-. Notice that the scales 
of the ordinates of the graphs differ within as well as between the two figures. 

In Fig. 9 we show the relative renormalization scale uncertainties of the singlet quark 
and gluon evolution, estimated using the smaller conventional range < /i 2 < 4/ij via 

max [q(x, (4 = \fj? f . . . - min [q(x, fi 2 r = |/j • ■ ■ 



Aq 



, q = H,g. (5.3) 



2 | average [q(x, $ = \n) . . . 4^)] | 

At large x this estimate yields uncertainties below 2% and 1% at NNLO for £ and g, 
respectively, representing improvements by more than a factor of three with respect to 
the corresponding NLO results. A clear improvement is also found for the singlet quark 
derivative £ at small x, persisting even down to x — 1CT 4 at fij = ~ 30 GeV 2 . For the 
gluon density, on the other hand, such a small- a; improvement will only occur at NNLO 
if p("p(x) is closer to the (small) approximation A than to the function B (cf. Fig. 4). 
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Figure 9: The relative /^-uncertainties of the NLO and NNLO predictions for £ and g 
at /z 2 = /i? ~ 30 GeV 2 , as estimated by the quantities AS and Ag defined in Eq. (5.3). 
The spikes close to x — 0.1 derive from zeros in the denominator of Eq. (5.3). 



We now turn to the perturbative expansion (2.16) for the singlet structure function 
F 2y s (henceforth simply denoted by F 2 ) and its Q 2 -derivative. The corresponding figures 
below refer to Q 2 = /i 2 ~ 30 GeV 2 and the input (5.1) and (5.2). For brevity the results 
will be displayed for /i r = [if = \i only. As above, the present NNLO uncertainties due to 
the incomplete information on the three-loop splitting functions are represented by the 
bands spanned by the NNLO^ and NNLO# approximations. 



The LO, NLO and NNLO results are compared in Fig. 10 for the standard scale 
choice /i 2 = Q 2 . Here the NLO and NNLO corrections to F 2 , shown in the left part of the 
figure, derive from the coefficient functions only. The positive corrections at large x stem 
from the non-singlet parts of the quark coefficient functions c^i, n = 1, 2 , which receive 
large contributions ~ [ln fc (l — x)/(l — x)] + , k = 0, . . . , 2n— 1, from soft-gluon emissions. 
At x = 0.8, for example, the NNLO term adds 15% to the NLO result, which in turn 
exceeds the LO expression ^S(x, Q 2 ) by 45%. The sizeable negative NNLO corrections 
at small x, on the other hand, are mainly due to the gluon contribution. In this region 
c 2^g{ x ) is dominated by the 1/x term arising from t-channel soft gluon exchange. Due 

' (2) 

to the convolution with the gluon density, however, the other contributions to c 2 g are as 



important, as already pointed out in refs. 0, |34j|. For our input (5.1) the total NNLO 
effect at 10~ 4 <x<10~ 2 amounts to about 40% of that without the (positive) 1/x terms. 
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Figure 10: Comparison of the LO, NLO and NNLO results for the singlet structure 
function F 2 (left part), its linear derivative F 2 ' = dF 2 /dlnQ 2 (right part, small x), and 
its logarithmic derivative F 2 = d In F 2 /d In Q 2 (right part, large x) at /j, 2 . = fi 2 = fi 2 = Q 2 . 
All results refer to the parton densities of Eq. (5.1) and a s = 0.2. 



Like the evolution of the singlet quark density, the Q 2 -derivative of F 2 shown in the 
right part of Fig. 10 is dominated by the quark contribution at large x, x > 0.3, and by the 
gluon contribution at small x, x < 0.03. It is thus convenient to consider the logarithmic 
derivative F 2 = dlnF 2 /dlnQ 2 in the former x-range, while in the latter region the linear 
derivative F 2 = dF 2 /d In Q 2 is more appropriate. Due to the effects of the quark coefficient 
functions discussed above, the LO and NNLO corrections to F 2 at large x are considerably 
larger than their counterparts for E illustrated in Fig. 5. The NNLO corrections rise from 
3% at x — 0.5 to 11% at x — 0.8, the corresponding NLO contributions amount to 18% 
and 24% of the LO results. The small positive gluon contribution to F 2 receives large 
corrections as well. At x — 0.5, for instance, it reaches 1.7%, 3.4% and 4.8% of the total 
|F 2 | at LO, NLO, and NNLO, respectively, in the latter approximation falling below 1% 
only at x ~ 0.7. F{ exhibits a NNLO effect of about 10% at 10~ 4 < x < 0.05, while the 
corresponding NLO/LO ratio rises from about 1.2 at x — 0.03 to 1.7 at x — 10~ 4 . This 
rather constant NNLO correction combines the effects of the coefficient functions (which 
dominate for x > 0.01, but decrease below) and the three-loop splitting functions (the 
impact of which rises towards very small x, cf. Fig. 5). The NNLO uncertainties due to 

(2) 

the incomplete information on P-j (x) are very similar to those for the quark evolution, 
i.e., they reach about ±2% at x — 10~ 3 and then sharply increase towards x^0. 
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The scale dependence of F 2 and of its Q 2 -derivatives is illustrated in Figs. 11 and 12 
for the same six x- values as in Figs. 7 and 8. Analogous to those figures /i is varied over 
the range | Q 2 < /i 2 < 8 Q 2 . The resulting estimates for the relative scale uncertainties 

max [fjx, » 2 = \Q 2 ... 4Q 2 )] - min [fjx, /i 2 = \Q 2 . . . 4Q 2 )] . 
J ~ 2|average[/(a;, y u 2 = ±Q 2 ...4Q 2 )]| ' 1 2 ' 2 ' 2 ' 

(5.4) 

are presented in Fig. 13. Except for the intermediate x-region 0.03 < x < 0.3 (where F 2 is 
quite stable already at NLO and F 2 is small) the NNLO results represent an improvement 
by a factor two or more at x > 10~ 4 for F 2 and 10~ 3 < x < 0.7 for F 2 \ 
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Figure 11: The dependence of the NLO and NNLO predictions for the singlet structure 
function F 2 (x, Q 2 = fj, 2 ) on the scale jj, for six representative values of x. 
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Figure 12: As Fig. 11, but for the logarithmic derivative F 2 (upper four x-bins) and the 
linear derivative F 2 ' (lower two values of x). 

It is also interesting to consider separately the dependence of F 2 and its derivatives on 
the renormalization scale ji r (keeping (j, 2 — Q 2 — fi 2 ) and on the factorization scale \if 
(keeping [i 2 = Q 2 = /J-y- ). For F 2 the resulting variations are comparable over the full 
arrange at NNLO. At x — 0.8, for example, they are both about half as large as the 
/^-dependence shown in Fig. 11. It is worth noting that the \i r - variation of F 2 at small-a; 
worsens by the transition from NLO to NNLO, despite the fact that a renormalization scale 
logarithm enters for the first time only at NNLO, see Eqs. (2.16) and (2.18). This is due 
to the large contribution of c^ g shown in Fig. 10. When both scales are varied, however, 
this deterioration is overcompensated by the clear improvement of the /i/-dependence. 
For the Q 2 -derivatives of Fig. 12 the effect of /i/ is much smaller than that of [i r , except 
for the NLO case at small x (where the large scale dependence of S of Fig. 9 is relevant). 
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Figure 13: The relative scale uncertainties of the NLO and NNLO predictions for F 2 (x, Q 2 ) 
and its linear (small x) and logarithmic (large x) Q 2 -derivatives at Q 2 ~ 30 GeV 2 , as 
estimated by the quantities AF 2 , AF 2 and AF 2 defined in Eq. (5.4). 



At NNLO the dependence on fif is relatively largest at very large x, but even at x = 0.8 
it does not exceed a quarter of the //-variation shown in Fig. 12. 

So far the NNLO effects have been discussed for the fixed initial parton densities 
E(x,/t 2 ) and g(x,fi 2 Q ) of Eq. (5.1) and the coupling constant a s (fi 2 ) = 0.2. We now 
finally address the impact of the NNLO terms on the determination of these perturbatively 
incalculable inputs from data on F 2 . For a simple estimate of the resulting difference 
between the NLO and the NNLO singlet parton densities and their respective theoretical 
uncertainties, we employ as 'data' the values of F 2 (x,Q 2 ) and F 2 (x,Q 2 ) for Q 2 = fi 2 ~ 
30 GeV 2 at the six x- values of Figs. 11 and 12, supplemented by x = 10 4 . I.e., the 
parameters of the input (5.1) (including a (1 + Ax) interpolation term for the gluon 
density, but imposing the momentum sum rule) are fitted at NNLO for /t 2 = Q 2 to 
the corresponding NLO results. The errors bands due to the scale dependence are then 
determined by repeating the NLO and NNLO fits to their respective fi 2 = Q 2 central 
values of F 2 and F 2 ' for fi 2 /Q 2 = 0.25, 0.5, 2 and 4. 

A critical point in analyses of F 2j s is the correlation between a s and the gluon density, 
which tends to strongly increase the //-variation of the fitted values of a s . Indeed, under 
the present conditions this variation is enhanced, in both NLO and NNLO, by about a 
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Figure 14: Comparison of the NLO and NNLO singlet quark and gluon densities and 
their scale-uncertainty bands as obtained by fitting F 2 and F 2 ' at Q 2 = fi 2 Q ~ 30 GeV 2 
for 10~ 4 < x < 0.8. The range of scales used for the error bands is 0.25 Q 2 < fi 2 < AQ 2 . 



factor of two with respect to an analogous analysis of the non-singlet structure function 
-p2,NS> thus preventing an accurate determination of a s from F 2t $ alone. In order to keep 
the scale variation of a s at a level realistic for analyses of data on electromagnetic DIS 
(where both -F 2 P ns anc ^ -^2 P s have been accurately measured), we therefore take over the 
values for a s (fi 2 ) from our previous non-singlet analysis ||: a s (fi 2 ) = 0.188, 0.200 and 
0.220, respectively, for /i 2 = 0.25 Q 2 , Q 2 and 4Q 2 and NLO. The corresponding NNLO 
values read 0.190, 0.194 and 0.202. It is clear from our discussion above that these scale 
uncertainties of a s are almost entirely due to the variation of the renormalization scale. 

The resulting NNLO central values and NLO and NNLO scale- uncertainty bands for 
fi 2 ) and g(x, fi 2 ) are presented in Fig. 14. The shift of Snnlo with respect to Enlo 
- an enhancement between 3% and 7% at small x, and a sizeable decrease at very large 
x reaching 13% at x — 0.8 - rather directly reflects the results for F 2 shown in Fig. 10. 
Likewise the reduction of the scale dependence to ±2% or less for 10 _3 < x < 0.5 closely 
follows the pattern of AF 2 in Fig. 13. The corresponding results for the gluon density, 
on the other hand, are considerably affected by the values of a s . In fact, the theoretical 
error band for g NNL o does not exceed 2% for 3-10~ 3 < x < 0.2. For x < 10~ 3 the present 
uncertainties of the three-loop splitting functions prevent an improvement on the rather 
accurate NLO determination of the gluon density. 
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6 Summary 



We have investigated the effect of the NNLO perturbative QCD corrections on the evo- 
lution of the unpolarized flavour-singlet quark and gluon densities, £(x,/z 2 ) and g(x,/j, 2 ), 
and on the most important singlet structure function, i^s^ar, Q 2 ). Our main new ingre- 
dients are approximate expressions for the rc-dependence P^ (x) of the singlet three- loop 
splitting functions, including quantitative estimates of their residual uncertainties. These 
parametrizations have been derived from the lowest four even-integer moments PA- (iV), 
N = 2, 4, 6 and 8, calculated in ref. ||, supplemented by the results of refs. ]TT|, [T2j for 



the leading small-x terms oc (l/s)lnx of Pffi, Pj$ and Consequently the differ- 

ences between our parametrizations illustrated in Figs. 3 and 4 oscillate at large x, as 
the difference between any two approximations obviously has four vanishing moments. 
At x < 10~ 2 the differences increase, because the [l/x)\nx terms unfortunately do not 
sufficiently dominate over the unknown less singular contributions like 1/x in the experi- 



mentally accessible region of small x, see also refs. [^4], [35| . 



Clearly our parametrizations represent only a temporary 'solution' of the problem of 
the three-loop splitting functions, which will be superseded sooner or later by an exact cal- 
culation. However, for two reasons the present approach is, over a wide region of x, much 
more effective than a brief look at Figs. 3 and 4 might suggest. First of all the splitting 
functions enter physical quantities only via convolutions with nonperturbative initial dis- 
tributions, which smoothen out the above-mentioned oscillations to a large extent. While 
this mechanisms leads to especially small uncertainties at x > 0.1, the convolutions extend 
the relevance of the large-x constraints || deep into the small- X region as exemplified in 
the discussion of Fig. 2. The second reason is that the series expansion for the evolution 
of the parton densities is very well converging — except, possibly, at x < 10~ 3 . Choosing 
for definiteness a s = 0.2 for the strong coupling constant (this corresponds to a scale 
between about 20 GeV 2 and 50 GeV 2 , depending on the precise value of a s (M|) ) the 
NNLO effect on the evolution of fj, 2 ) and g(x,/i 2 ) amounts to less than 2% and 1%, 
respectively, at x > 0.2. Hence the net uncertainties due to the incomplete information 
on the Pj-j (x) are absolutely negligible in this region as shown in Figs. 5 and 6. Of course, 
these uncertainties rise towards small x, but they exceed ±2% only below x ~ 10~ 3 (or 
a few times this number, if much lower scales are involved), leaving a comfortably large 
region of safe applicability of our results to processes at the Tevatron and the LHC. 

At x > 0.2 the NNLO renormalization-scale variation (for the conventional interval 
\l j2 f < A*r — amounts to less than 2% and 1%, respectively, for dT./d In fx 2 and 

dg/d In fj, 2 at our reference point n 2 ~ 30 GeV 2 . The corresponding numbers at x = 10~ 3 
read 5% and 3%. Except for the gluon evolution at this latter value of x, these results 
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represent an improvement on the NLO evolution by a factor of three or more. Taking 
into account also the rapid convergence at \i T — {if and our corresponding findings for the 
non-singlet sector in ref. |J, we expect that terms beyond NNLO will affect the parton 
evolution by less than 1% at large x and less than 2% down to x ~ 10~ 3 for a s ~ 0.2. 

Due to the additional effect of the two- loop coefficient functions [Q, [|, the structure 
function F 2 ^s and its scaling violations receive considerably larger NNLO corrections at 
x > 10~ 2 . In fact, keeping only the coefficient functions (and omitting the three- loop 
splitting functions Pj? (x) altogether) forms quite a good approximation in this x-range. 
As shown in Fig. 10, the NNLO corrections for both F 2 and its Q 2 -derivative are partic- 
ularly large at very large x, e.g., 15% for F 2y s and 11% for d In F 2y g/d In Q 2 at x — 0.8 for 
our reference scale Q 2 ~ 30 GeV 2 . This is an effect of the large soft-gluon terms in the 
quark coefficient functions which are not special to the singlet case considered here. While 
the gluon contributions dominates the sizeable (up to about 7%) NNLO corrections at 
small x, their effect is suppressed at large x, especially for the absolute size of F 2 ^ S - It is 
worth noting, however, that that the gluon contribution to dF 2 s/d\aQ 2 at x = 0.5 still 
amounts to 5% at NNLO (40% more than at NLO), an effect large enough to jeopardize 
analyses which apply a purely non-singlet formalism to the data on the proton structure 
function F 2 (x, Q 2 ) in the region x > 0.3. The accuracies of both F 2 and its scaling viola- 
tions, as estimated by the scale dependence, are considerably improved (by more than a 
factor of two over a wide region) by the inclusion of the NNLO terms. Consequently the 
same applies to the theoretical accuracy of determinations of the singlet quark and gluon 
densities from data on F 2t s and dF 2 s/dlnQ 2 at Q 2 ~ 30 GeV 2 illustrated in Fig. 14: 
uncertainties of less than 2% from the truncation of the perturbation series are obtained 
for the quark density at 10 _3 <a;<0.5 and for for the gluon density at 3 • 10 _3 <x<0.2. 

Our results of ref. || and the present paper complete, if only approximately at x > 10~ 3 , 
the theoretical prerequisites for NNLO analyses of structure functions in DIS and total 
cross sections of Drell-Yan processes. Further progress at large x — especially a fur- 
ther improvement on the theoretical accuracy of NNLO a s -determinations from structure 
functions |36|, |37|, || — can be obtained by including the 0(a 3 ) coefficient functions, par- 



ticularly in the non-singlet sector dominating the extraction of a s . In fact, results on 
these functions are available from the fixed-moment calculations of refs. [0, || and from 
soft-gluon resummation p8[ . We will address this issue in a forthcoming publication. 
Progress towards the important HERA small-x region of x < 10~ 3 at moderate to low Q 2 , 
however, definitely requires the full calculation of the three-loop splitting functions. 

(2) (2) 

Fortran subroutines of our parametrizations of c a j (x) (a = 2, L; i = q, g) and (x) 
can be obtained via email to neerven@lorentz.leidenuniv.nl or avogt@lorentz.leidenuniv.nl. 
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Appendix: The singlet coefficient functions in AT-space 

The Mellin transforms (2.12) of the parametrizations (3.3)-(3.6) of the two-loop singlet 
coefficient functions for Fl and F 2 are given in terms of the integer- iV sums Si(N) and 
their complex- iV analytic continuations 

N 1 C-1V 
Si = Sm = E ] j=ai)-j r ^y^- 1) (N+l) . (A.l) 

Here ((1) represents the Euler-Mascheroni constant, and ((l> 1) Riemann's (^-function. 
The lih logarithmic derivative of the T-function can be readily evaluated using the 

asymptotic expansion for Re TV > 10 together with the functional equation. 

The moments (2.12) of the pure singlet coefficient functions (3.3) and (3.4) are given by 



( 2 ) , «. Ar \t 15.94 37.092 26.364 5.212 x 

c Lnv = Nf<[ + + )s 1 

L ' PSK ' f \\ N N + l N + 2 N + 3J 

2.370 0.842 28.09 26.38 3.040 65.182 

+ 77^ 7777" + —TT- + TT7 77^ + 



N-l N 3 N 2 N (N + l) 3 (iV+1) 2 
88.678 26.364 91.756 5.212 27.0881 
N +1~ (N + 2) 2 + N + 2 + (N + 3) 2 ~ N + 3) ^ 



and 



(2) , Ar . Ar f/49.702 27.802X /49.5 30.23 27.903 ... 
c?h(N) = N f \(— —--—)S 3 +{ — + — -——, S>. 



/ 1 1 \ / „ „ /49.5 30.23 28.206 \ „ 



+ 



(N + l) 3 , 

5.290 25.86 4.172 121.205 114.519 83.406 



N — 1 iV 4 N 3 N 2 N (N+ iy 



45.4003 33.1769) /A , 



27 



The corresponding expressions for the gluon coefficient functions (3.5) and (3.6) read 
(2)/ArX Ar f / 94.74 1017.06 49.20 \ „ Z94.74 143.94 49.20 \ o2 

864.8 873.12 963.2 98.40 98.40 \ 

~W~ + (N + l) 2 + JV + 1 + (N + 2) 2 ~ iV"+2 J 1 

5.333 39.66 5.333 2154.24 1002.86 1909.768 

+ 



N-l N 2 N (iV + 1) 3 (iV + 1) 2 N+l 
98.40 98.40 ) /k t . 

+ titt^t, ( A -4) 



(N + 2) 3 (iV + 2) 5 



and 



4>) = N f \{— — + — )s 3 + ( — -- + — — 2 )s 2 




628.2 

(ivTT)' 

628.2 
(iV+1) 2 



S 2 (A.5) 



s , _ - 215.845 _ 209^4 \ 3 
iV 3 N {N+iyJ 1 V AT N+l) 1 

31.914 118.96 2097.4 4938.34 1256.4 _ 271 ] 



JV-1 iV 4 N 3 N 2 N (N + iy 



Let us finally briefly elaborate on the derivation of Eq. (4.7) from Eq. (4.5). The exact 
expression §, g| for the residue of the N = 1 pole of t^j(iV) reads iV> [344/9 - 16 ((2)]. 
Recalling that the Mellin transform of (1/x) \nx is given by — 1/(N— l) 2 , this information 
combined with the corresponding term 16/[3(iV— 1)] of the LO splitting function is 
sufficient to obtain Eq. (4.7) from Eqs. (4.5) and (4.6). The latter relation simply arises 
since the other contributions to the scheme transformation (2.24) do not lead to any 
l/(iV-l) 2 terms. 
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